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Ecosystem  Predictions  with  Approximate  vs.  Exact  Light  Fields 


1.  Executive  Summary 

Coupled  physical-biological-optical  ecosystem  models  often  use  sophisticated  treatments  of  their 
physical  and  biological  components,  while  oversimplifying  the  optical  component  to  the  detriment 
of  the  ecosystem  predictions.  To  bring  optical  computations  up  to  the  standard  required  by  recent 
ecosystem  models,  I  previously  developed  (as  a  part  of  the  Hyperspectral  Coastal  Ocean  Dynamics 
Experiment  (HyCODE)  under  ONR  contract  N000 1 4-05 -M-0 146)  a  computationally  fast  numerical 
model  named  EcoLight,  derived  from  the  HydroLight  numerical  model,  which  computes  accurate 
in-water  irradiances  for  use  in  biological  primary  production  calculations.  In  the  present  work,  my 
colleagues  and  I  incorporated  EcoLight  into  a  coupled  physical-biological  model.  We  then  compared 
predicted  chlorophyll  concentrations  and  nutrients  for  an  idealized  open-ocean  simulation  of  mid¬ 
latitude  Case  1  waters  for  approximate  analytical  vs.  accurate  numerical  (EcoLight)  irradiance 
calculations. 

In  these  simulations,  the  predicted  chlorophyll  concentrations  differed  by  a  few  tens  of  percent 
over  the  course  of  10-year  simulations  owing  to  the  differences  in  the  approximate  vs.  accurate 
irradiances.  This  speaks  well  of  the  approximations  used  within  the  analytic  light  model  for  Case 
1  water.  More  importantly,  the  EcoLight  code  runs  with  no  more  than  a  30%  increase  in  total  run 
time,  while  providing  several  advantages  than  cannot  be  obtained  with  an  analytic  light  model. 
These  advantages  are 

•  EcoLight  does  a  much  better  job  of  computing  irradiances  near  400  and  beyond  00  nm, 
and  at  great  depths.  These  differences  can  have  significant  effects  on  how  particular 
phytoplankton  functional  groups  evolve  in  time  and  at  depth. 

•  EcoLight  can  model  any  water  body,  including  Case  2  water  and  optically  shallow  waters 

for  which  bottom  reflectance  can  substantially  increase  the  irradiance  available  for 
photosynthesis  and  water  heating. 

•  EcoLight  provides  the  outputs  necessary  for  ecosystem  validation  using,  for  example, 
satellite-derived  water- leaving  radiances  or  remote- sensing  reflectances,  without  the  necessity 
of  using  a  chlorophyll  algorithm  to  convert  satellite  data  into  an  ecosystem  variable. 

•  EcoLight  provides  a  means  for  improved  ecosystem  model  initialization  and  data 
assimilation  directly  in  terms  of  the  available  data  (IOPs,  irradiances,  remote  sensing 
reflectances),  again  without  an  intermediate  step  involving  a  data-to-chlorophyll  conversion 
algorithm. 

•  EcoLight  costs  at  most  30%  extra  in  the  total  simulation  run  time.  This  is  a  small 
computational  price  to  pay  for  its  advantages. 
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2.  Background 


The  fundamental  measure  of  light  energy  in  an  aquatic  system  is  the  spectral  radiance,  which  in 
horizontally  homogeneous  waterbodies  is  a  function  of  depth,  direction,  and  wavelength.  However, 
for  ecosystem  modeling  the  directional  information  contained  in  the  radiance  is  usually  irrelevant, 
because  water  constituents  such  as  phytoplankton  and  dissolved  substances  are  assumed  equally 
likely  to  interact  with  light  regardless  of  its  direction  of  travel  or  state  of  polarization.  Therefore, 
the  spectral  scalar  irradiance  as  a  function  of  depth  z  and  wavelength  e,  EJz.c),  is  the  fundamental 
radiometric  quantity  necessary  for  predictions  of  aquatic  primary  productivity,  heating  of  water,  and 
photochemical  reactions.  When  modeling  photosynthesis,  which  depends  on  the  number  of  photons 
absorbed,  it  is  necessary  to  multiply  Ea(z,e )  by  cl tu\  where  h  is  Planck's  constant  and  c  is  the  speed 
of  light.  This  converts  the  energy  units  of  the  spectral  irradiance,  W  m 2  nm  ',  to  quantum  units, 
photons  s'1  m  2nm 

A  wavelength-integrated  measure  of  the  total  light  available  for  photosynthesis,  the 
Photosynthetically  Available  Radiation  (PAR),  has  historically  been  used  in  simple  ecosystem 
models.  Simple  analytical  models  are  commonly  employed  for  estimating  the  dependence  of  PAR 
on  time  (diurnal  to  seasonal  changes),  depth,  and  water  properties.  These  models  are  often 
parameterized  by  the  chlorophyll  concentration.  Climatological  data  are  often  used  to  provide  the 
time  dependence  of  the  (down  welling)  surface  irradiance.  Traditionally,  coupled  ecosystem  models 
have  incorporated  simple  analytical  formulas  to  predict  PAR(z )  from  a  given  chlorophyll  profile  and 
from  an  estimate  of  PAR  at  the  sea  surface.  Although  these  analytical  PAR  models  are 
computationally  fast,  they  can  produce  estimates  that  differ  by  a  factor  of  three  near  the  sea  surface 
and  by  a  factor  of  ten  at  the  bottom  of  the  euphotic  zone  (Zielinski  et  al.,  1998).  Indeed,  the  depth 
of  the  euphotic  zone,  if  defined  as  the  depth  where  PAR  decreases  to  one  percent  of  its  surface  value, 
differs  by  almost  a  factor  of  two  among  these  models.  Errors  of  this  magnitude  are  unacceptable  for 
quantitative  predictions  of  primary  productivity  or  upper-ocean  thermal  structure.  Additional 
inaccuracies  can  be  found  in  the  application  of  the  PAR  formulas,  such  as  the  use  of  downwelling 
plane  irradiance  Ed  as  a  proxy  for  the  scalar  irradiance  E0 . 

Although  simple  PAR  models  are  computationally  convenient,  any  ecosystem  model  based  on 
PAR  is  limited  in  its  realism,  regardless  of  how  accurately  PAR  is  computed.  Different  species  of 
phytoplankton,  or  even  the  same  species  under  different  environmental  conditions,  have  different 
pigment  suites  and  thus  respond  differently  to  the  same  light  field.  Therefore,  any  biological  model 
that  attempts  to  describe  competition  between  different  species  or  functional  groups  having  different 
pigments  must  use  spectral  irradiance,  not  a  wavelength-integrated  measure  of  the  light. 

The  Ecosystem  Simulation  (EcoSim)  model  of  Bissett  et  al.  (1999a,  1999b)  includes  four 
phytoplankton  functional  groups,  each  with  a  characteristic  pigment  suite  that  changes  over  time 
with  light  and  nutrient  conditions.  To  model  light  effects  on  changes  within  and  competition 
between  these  groups,  EcoSim  uses  an  approximate  analytical  model  for  the  spectral  scalar 
irradiance  in  Case  1  waters  (described  below).  Although  the  EcoSim  spectral  irradiance  model  is 
sufficiently  fast  for  use  in  coupled  physical-biological  ecosystem  simulations  involving  many  grid 
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points  and  time  steps,  it  does  rest  on  several  approximations  that  limit  its  accuracy  even  in  Case  1 
waters,  and  it  is  not  applicable  to  optically  shallow  waters. 

There  is  today  much  interest  in  the  modeling  of  coastal  waters,  which  are  often  Case  2  due  to 
resuspended  sediments  or  terrigenous  particles  and  dissolved  substances,  and  which  may  be  shallow 
enough  for  bottom  reflectance  to  make  a  significant  contribution  to  the  scalar  irradiance.  If  EcoSim 
or  similar  models  (e.g.,  Fujii,  et  al.,  2007)  are  to  accurately  simulate  such  optically  complex  waters, 
they  must  incorporate  spectral  irradiance  models  that  are  computationally  fast,  accurate,  and 
applicable  to  any  water  body  or  bottom  condition.  To  address  this  need,  I  developed  a  specialized 
version  of  the  HydroLight  radiative  transfer  model,  called  EcoLight,  which  is  designed  to  bring  the 
optical  component  of  ecosystem  models  up  to  the  level  of  sophistication  found  in  the  latest  physical 
and  biological  models. 


3.  Work  Completed 

The  basic  EcoLight  code  was  developed  with  previous  ONR  funding  and  preliminary  ecosystem 
simulations  were  made  to  show  the  importance  and  computational  feasibility  of  accurate  irradiance 
calculations.  The  current  work  incorporated  EcoLight  into  a  coupled  physical-biological  model  and 
compared  its  performance  with  the  analytical  spectral  irradiance  model  used  in  EcoSim.  I  first 
briefly  describe  the  physical  and  biological  models,  and  then  EcoLight  and  the  coupled  model.  I  then 
compare  the  predictions  made  for  an  idealized  ecosystem  simulation  when  using  the  original  EcoSim 
light  model  with  those  made  using  various  EcoLight  options. 

The  work  reported  here  was  performed  in  collaboration  with  Lydia  Sundman  of  Sundman 
Consulting,  Paul  Bissett  of  the  Florida  Environmental  Research  Institute,  and  Bronwyn  Cahill  of 
Rutgers  University,  who  were  supported  on  this  contract.  All  have  contributed  to  this  final  report. 
This  report  is  condensed  from  a  draft  journal  article  (Mobley  et  al.,  2009),  which  is  in  preparation 
for  submission  to  Biogeosciences  Discussions. 

3.1  The  ROMS  physical  model 

The  physical  model  used  here  is  a  version  of  the  Regional  Ocean  Modeling  System  (ROMS) 
described  by  Shchepetkin  and  McWilliams  (2005;  see  also  Fringer  et  al.,  2006  and 
http://marine.mtgers.edu/po/index.php?model=roms).  The  ROMS  model  is  a  curvilinear-coordinate 
(terrain-following),  free-surface,  primitive  equation  model  designed  for  prediction  of  physical 
oceanography  quantities  in  coastal  waters.  The  present  study  used  a  6x6  horizontal  grid  version  of 
ROMS  3.0  with  periodic  lateral  boundary  conditions.  This  spatially  limited  ROMS  version  was 
chosen  to  minimize  run  times  during  the  EcoLight  code  development  and  verification  simulations. 
For  the  simulations  presented  below,  the  grid  is  centered  off  the  continental  shelf  of  the  eastern 
United  States  near  72.8  deg  West  and  38.8  deg  North.  The  horizontal  grid  resolution  is 
approximately  13  km.  The  vertical  grid  has  30  points  covering  the  upper  202  m  of  the  water  column; 
the  vertical  resolution  ranges  from  approximately  2  m  near  the  sea  surface  to  15  m  at  depth. 
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3.2  The  EcoSim  biological  model 

The  Ecosystem  Simulation  (EcoSim)  biological  model  (Bissett  et  al.,  1999a,  1999b)  was 
developed  for  simulations  of  carbon  cycling  and  biological  productivity.  This  model  includes  four 
phytoplankton  functional  groups,  defined  according  to  their  pigment  suites  (small  diatoms,  large 
diatoms,  dinoflagellates,  and  synechococcus).  Each  functional  group  has  a  unique  set  of  accessory 
pigments,  which  varies  with  the  group  carbon-to-chlorophyll  a  ratio,  C:Chl.d.  Pigment  packaging  and 
accessory  pigment  concentration  are  functions  of  the  chlorophyll  a  concentration  within  each 
functional  group.  The  chlorophyll  a  content  and  other  properties  of  each  functional  group  evolve 
with  the  light  history  and  nutrient  status  of  the  group.  The  model  also  includes  components 
describing  dissolved  and  particulate  organic  and  matter,  bacteria,  and  detritus.  The  interactions 
between  these  components  describe  autotrophic  growth  of  and  competition  between  the  four 
phytoplankton  groups,  differential  (non-Redfield  ratio)  carbon  and  nitrogen  cycling,  nitrification, 
grazing,  and  air-sea  exchange  of  C02.  The  initial  application  of  EcoSim  to  predictions  of  seasonal 
cycles  of  carbon  cycling  and  phytoplankton  dynamics  in  the  Sargasso  Sea  showed  that  its  predictions 
were  consistent  with  measurements  of  various  biological  and  chemical  quantities  at  the  Bermuda 
Atlantic  Time-series  Study  station  (Bissett  et  al.,  1999a). 

An  important  component  that  was  missing  from  the  original  EcoSim  code  was  nutrient  recycling. 
For  the  present  simulations,  subroutines  were  added  to  allow  for  vertical  particle  flux  and  restoration 
of  nutrient  fluxes  that  “fall  through”  the  computational  bottom  back  into  the  nutrient  pool.  Fecal 
material  and  phytoplankton  now  sink  in  EcoSim.  When  these  hit  the  bottom,  the  flux  of  material 
from  each  component  out  of  the  system  is  converted  directly  into  the  appropriate  nutrient  and 
restored  back  into  the  nutrient  pool.  The  result  of  this  effort  is  that  we  can  now  obtain  a  nearly  stable 
annual  cycle  over  multi-year  simulation  periods. 

The  absorption  spectra  of  the  phytoplankton  functional  groups  change  with  light  and  nutrient 
adaptation.  The  four  groups  therefore  respond  differently  to  various  wavelengths  of  the  available 
light,  and  each  particular  group  responds  differently  over  time.  EcoSim  requires  spectral  irradiances 
at  5  nm  bandwidths  between  400  and  700  nm  in  order  to  model  the  changes  within  each  functional 
group  and  competition  between  them.  The  use  of  spectral  irradiance  rather  than  broadband  PAR  in 
modeling  phytoplankton  dynamics  is  a  distinguishing  feature  of  EcoSim. 

EcoSim  uses  the  RADTRAN  atmospheric  model  (Gregg  and  Carder,  1990)  to  compute  spectral 
downwelling  plane  irradiances  just  beneath  the  sea  surface,  Ed( 0,e).  These  irradiances  can  be 
rescaled  to  match  given  irradiance  values  if  the  model  is  being  driven  by  measured  or  climatological 
sky  conditions.  These  subsurface  spectral  downwelling  plane  irradiances  are  then  propagated  to 
depth  using 


E^zX)  =  EJfiX)  exp 


fKd(z/,X)dZ/ 

0 


(i) 
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and  a  simple  model  for  Kd, 


a(z,X)  +  bb(z,X) 

- -  /  vx - •  (2> 

Here  a( z,e)  is  the  total  absorption  coefficient  (the  sum  of  absorption  by  pure  water  and  the  various 
particulate  and  dissolved  components),  and  bh{ z,e)  is  the  total  backscatter  coefficient.  The 
phytoplankton  absorption  is  obtained  from  the  concentrations  of  the  functional  groups  and  their 
chlorophyll-specific  absorption  spectra.  The  backscatter  and  total  scatter  coefficients  are  obtained 
from  chlorophyll-dependent  models  for  Case  1  waters,  using  the  total  chlorophyll-a  concentration. 
The  mean  cosine  for  downwelling  irradiance,  i  d(z,e),  is  itself  modeled  by  a  simple  function  that 
merges  estimates  of  the  near-surface  and  asymptotic-depth  mean  cosines  (Bissett  et  al.,  1999b,  Eqs. 
18-22).  Finally,  the  needed  scalar  irradiance  Ea( z,e)  is  obtained  from  the  computed  Ed( z,e)  via  the 
approximation  E0( z,e)  =  Ed(z,'6)Kd(z,e)/a(z,e). 

The  biology  is  updated  at  each  time  step  and  grid  point  using  the  analytic  formulas  for  the  scalar 
irradiance  just  described.  However,  the  irradiances  computed  within  EcoSim  do  not  feed  back  to 
the  ROMS  code  which,  for  programming  simplicity  when  merging  the  codes,  retains  its  original 
short  and  long-wave  light  parameterization  for  mixed-layer  heating  calculations.  Thus  the  physical 
model  influences  the  biology  via  temperature  and  mixing,  but  the  optical  model  employed  within 
EcoSim  does  not  influence  the  physical  model.  This  simplification  was  made  in  the  present  study 
to  avoid  alterations  to  the  ROMS  code. 

3.3  The  Ecolight  irradiance  model 

Numerical  models  that  solve  the  exact  radiative  transfer  equation  (RTE)  are  currently  available 
(Mobley,  et  al.,  1993)  and  can  provide  highly  accurate  predictions  of  the  light  field,  in  particular  the 
scalar  irradiance  as  a  function  of  depth  and  wavelength.  However,  the  computational  times  required 
by  these  exact  numerical  models  are  far  too  great  for  their  inclusion  in  large  coupled  models  that 
need  irradiance  predictions  at  many  spatial  locations  and  times  of  day  during  multi-year  ecosystem 
simulations. 

The  HydroLight  radiative  transfer  model  (Mobley,  et  al.,  1993;  Mobley,  1994; 
www.hydrolight.info)  provides  an  accurate  solution  of  the  unpolarized  RTE  for  any  water  body, 
given  the  inherent  optical  properties  (IOPs,  namely  the  absorption  and  scattering  properties)  of  the 
water  body,  the  incident  sky  radiance,  and  the  bottom  reflectance  (in  finite-depth  waters).  Although 
the  standard  version  of  Hydrolight  is  computationally  very  efficient  for  predicting  full  radiance 
distributions,  its  run  times  are  still  much  too  long  for  use  in  ecosystem  models. 

I  therefore  developed  a  highly  optimized  version  of  Hydrolight  4.2,  called  EcoLight,  which 
computes  irradiances  as  a  function  of  depth  and  wavelength  for  any  water  body,  with  the  same 
accuracy  as  HydroLight.  Related  quantities  such  as  the  irradiance  reflectance,  nadir- viewing  remote¬ 
sensing  reflectance,  diffuse  attenuation  functions,  and  mean  cosines  are  also  computed  by  EcoLight. 
Although  those  ancillary  quantities  are  not  needed  as  inputs  to  most  ecosystem  models,  they  can  be 
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of  use  in  ecosystem  model  validation  and  data  assimilation.  Quantities  such  as  the  remote-sensing 
reflectance  are  not  available  from  any  of  the  simple  analytic  light  models. 

Ecosystem  models  require  only  the  scalar  irradiance  as  a  function  of  depth  and  wavelength, 
which  is  computed  from  an  azimuthal  integration  of  the  radiance.  This  means  that  the  RTE  can  be 
azimuthally  averaged  to  obtain  an  equation  for  the  azimuthally  averaged  radiance  as  a  function  of 
polar  angle.  The  numerical  solution  of  the  resulting  RTE  is  much  faster  than  in  HydroLight. 
Various  other  simplifications  to  the  EcoLight  code  were  made.  The  most  important  are  the  use  of 
a  piecewise  homogeneous  water  column,  solution  of  the  RTE  to  different  depths  at  different 
wavelengths,  and  wavelength  skipping. 

The  ROMS-EcoSim  code  models  the  water  column  as  a  stack  of  homogeneous  layers  of  variable 
thickness.  Therefore  the  IOPs  within  a  given  water  layer  are  independent  of  depth.  The  EcoLight 
code  takes  advantage  of  the  depth  independence  of  the  IOPs  within  a  layer  to  reduce  the  computation 
needed  to  solve  the  RTE  for  the  depth  dependence  of  the  irradiances  within  each  layer. 

HydroLight  solves  the  RTE  to  a  user-specified  geometric  depth,  which  is  the  same  for  every 
wavelength  and  must  be  chosen  before  the  run  is  started.  This  becomes  computationally  expensive 
at  wavelengths  where  the  IOPs  are  large  (e.g.,  due  to  water  absorption  at  red  wavelengths), 
corresponding  to  a  large  optical  depth  for  a  given  geometric  depth.  EcoSim  requires  the  spectral 
scalar  irradiance  only  to  the  bottom  of  the  euphotic  zone,  below  which  it  does  not  perform  primary 
production  calculations.  Therefore,  it  is  necessary  to  compute  the  irradiance  only  to  that  depth, 
which  varies  with  the  biological  state  of  the  ocean  and  cannot  be  predetermined. 

In  EcoLight  the  goal  is  to  solve  the  RTE  to  the  shallowest  depth  possible  and  then  to  extrapolate 
the  scalar  irradiance  to  greater  depths.  To  determine  the  depth  to  which  the  RTE  is  to  be  solved  at 
a  particular  wavelength,  an  estimate  is  first  made  of  the  depth  za  to  which  the  scalar  irradiance  at  the 
current  wavelength  will  decrease  to  a  factor  F0  of  the  surface  value  (e.g.,  F0  =  0.1,  corresponding  to 
the  10%  light  level).  The  RTE  is  then  solved  exactly  to  depth  z,0.  The  computed  scalar  irradiance 
Ea( zQ,e)  is  then  objectively  extrapolated  to  all  depths  below  z0  using  a  variant  of  Eq.  (1),  namely 

z 

Eq(z,X)  =  E0(z0,X)  exp  -  J  a(z  ^ ,X)  dz  ^ 


HydroLight  numerical  simulations  show  that  the  scalar  irradiance  varies  with  depth  in  correlation 
with  the  depth  variation  of  the  total  absorption  coefficient  a(z',e),  once  zD  is  below  the  depth  at  which 
surface  boundary  effects  are  important.  This  is  to  be  expected  because  Kd  in  Eq.  (1)  is  a  very 
“absorption-like”  parameter,  and  Ka  is  close  to  Kd  away  from  the  sea  surface  (becoming  equal  to  Kd 
at  great  depths). 

There  is  also  an  EcoLight  option  to  solve  the  RTE  at  only  some  wavelengths,  and  to  obtain  the 
irradiances  at  the  unsolved  wavelengths  by  interpolation.  For  example,  EcoLight  can  solve  the  RTE 
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at  every  nth  EcoSim  wavelength  in  =  2,  3,  ...)  and  then  estimate  the  irradiances  at  the  unsolved 
wavelengths  by  interpolation  between  the  computed  wavelengths.  Solving  the  RTE  only  at  every 
nth  wavelength  gives  a  factor  of  n  decrease  in  the  EcoLight  run  time,  all  else  being  equal. 

In  summary,  EcoLight  takes  the  following  philosophy.  It  is  necessary  to  solve  the  RTE  in  order 
to  incorporate  the  effects  of  the  surface  boundary  conditions  and  to  account  for  all  IOP  effects. 
However,  once  an  accurate  value  of  E0( za  ,e)  has  been  computed  to  some  depth  za  deep  enough  to 
be  free  of  surface  boundary  effects,  it  is  not  necessary  to  continue  solving  the  RTE  to  greater  depths, 
which  is  computationally  expensive.  As  shown  in  Section  4,  it  is  possible  to  extrapolate  the 
accurately  computed  upper-water-column  irradiances  to  greater  depths  and  still  obtain  irradiances 
that  are  acceptably  accurate  for  ecosystem  predictions.  Likewise,  as  shown  below,  it  is  not  necessary 
to  solve  the  RTE  at  every  EcoSim  wavelength  in  order  to  obtain  acceptably  accurate  irradiances  at 
the  needed  wavelength  resolution. 

3.4  The  ROMS-EcoSim-EeoLight  coupled  model 

The  ROMS  and  EcoSim  models  were  previously  coupled  by  P.  Bissett,  and  that  code  served  as 
the  starting  point  for  the  nutrient  recycling  modifications  and  the  incorporation  of  EcoLight.  The 
analytic  irradiance  model  used  in  the  standard  EcoSim  code  was  replaced  by  a  call  to  the  EcoLight 
subroutine.  EcoSim  passes  EcoLight  the  current  total  IOPs  as  functions  of  depth  and  wavelength, 
atmospheric  conditions  (as  needed  by  RADTRAN,  which  is  also  used  by  EcoLight  to  compute  the 
spectral  irradiance  incident  onto  the  sea  surface),  wind  speed,  grid  depths,  and  other  information 
needed  by  EcoLight.  After  solving  the  RTE,  EcoLight  returns  the  scalar  irradiance  E0(z„e)  to  EcoSim 
for  use  in  its  primary  production  calculations.  The  diffuse  attenuation  function  Kd,  mean  cosine  i  d, 
remote-sensing  reflectance,  downwelling  plane  irradiance  Ed,  and  various  other  quantities  are  also 
returned  to  EcoSim.  Although  these  quantities  are  not  needed  for  EcoSim’ s  biological  computations 
when  EcoLight  is  used,  they  are  of  interest  for  comparison  with  the  analytic  values  and  for  model 
validation  when  comparing  ecosystem  predictions  with  actual  measurements.  These  quantities  are 
saved  to  files  for  post-run  analysis.  In  the  current  code,  the  EcoLight-computed  irradiances  are  used 
only  within  EcoSim;  ROMS  still  uses  its  original  irradiance  model  for  its  thermodynamics 
calculations.  Although  the  ROMS  parameterization  of  solar  transmission  into  the  water  column  is 
oversimplified  (Paulson  and  Simpson,  1977;  Cahill  et  al.,  2008),  it  was  retained  here  to  avoid 
modifications  to  the  ROMS  code  itself. 
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4.  Results 


4.1  Simulations 

The  ROMS-EcoSim  and  ROMS-EcoSim-EcoLight  models  were  run  in  various  modes  for  ten- 
year  simulations  to  determine  the  differences  in  ecosystem  evolution  owing  to  differences  in  the 
irradiance  calculations,  all  else  (ecosystem  structure,  vertical  mixing,  external  forcing)  being  held 
the  same. 

As  noted  in  §3.1,  the  6x6  ROMS  grid  used  here  was  located  off  the  east  coast  of  the  U.S.  The 
corresponding  external  forcing  for  solar  irradiance  and  surface  wind  speed  as  used  in  ROMS  was 
obtained  from  a  one-year  time  series  of  measured  values  from  ECMWF  40-year  reanalysis  of  gridded 
atmospheric  data  at  six-hour  resolution.  We  used  both  daily  averages  and  hourly  values  derived  from 
these  data  in  our  simulations.  The  RADTRAN  clear-sky  irradiance  computations  used  within 
EcoSim  were  driven  the  latitude,  longitude,  and  time  of  day,  and  the  RADTRAN  clear-sky  spectral 
irradiances  were  rescaled  to  match  the  measured  wavelength-integrated  values  used  by  ROMS .  The 
ROMS-computed  circulation  and  mixing  are  identical  in  all  runs  because  the  irradiances  computed 
and  used  within  EcoSim  for  its  biological  calculations  are  not  passed  back  to  ROMS  for  its  water¬ 
heating  calculations.  The  differences  in  the  predicted  ecosystem  evolution  are  thus  due  to 
differences  in  the  irradiances  computed  within  EcoSim,  either  by  its  default  analytic  model  or  by 
EcoLight,  and  used  there  for  its  primary  production  calculations. 

Given  the  geographical  location  of  the  ROMS  grid,  its  periodic  lateral  boundary  conditions,  and 
the  external  forcing,  the  present  simulations  can  be  thought  of  as  modeling  idealized  mid-latitude, 
open-ocean,  Case  1  water.  We  did  not  attempt  here  to  model  any  specific  ecosystem,  which  would 
require  a  fully  3D  spatial  grid  tailored  to  the  boundary  conditions  of  the  location  being  modeled. 

This  annual  cycle  of  external  forcing  was  repeated  for  each  year  of  the  simulation.  The  long-term 
annual  cycle  of  the  ecosystem  variables  should  thus  depend  only  on  the  ecosystem  structure,  i.e.,  on 
the  form  of  the  interactions  between  the  biological  functional  components,  on  the  nature  of  the 
vertical  mixing  used  in  the  physical  model,  and  on  the  values  of  parameters  that  describe  strengths 
of  various  interactions.  The  assumed  ecosystem  initial  conditions  and  numerical  transients  should 
affect  only  some  initial  time  period  of  a  simulation,  after  which  the  ecosystem  should  reproduce  an 
annual  cycle  that  is  independent  of  the  initial  conditions. 

4.2  Model  initial  conditions 

The  coupled  ROMS-EcoSim  model  with  its  analytical  light  calculations  was  first  run  to  examine 
its  behavior  over  long  simulation  times.  To  test  the  independence  of  the  final  ecosystem  state  from 
the  initial  conditions,  this  model  was  run  with  baseline  initial  (at  time  0)  concentration  vs.  depth 
profiles  for  the  biological  functional  groups.  The  baseline  initial  concentrations  for  the  four 
phytoplankton  functional  groups  are  shown  in  Fig.  1.  Those  concentrations  were  then  halved  and 
doubled  to  obtain  other  sets  of  initial  conditions.  Figure  2  shows  the  results  for  the  total  chlorophyll 
concentration  at  a  depth  of  1  m  for  the  first  year  of  the  simulations.  Halving  and  doubling  the  initial 


concentrations  affects  the  predicted  total  chlorophyll  values  for  only  the  first  few  months  of  a 
simulation.  Larger  perturbations  in  the  initial  conditions  require  longer  to  converge,  but  the 
predictions  do  approach  the  same  annual  cycle  in  later  years.  The  ecosystem  final  state  is  thus 
independent  of  the  initial  conditions  of  a  simulation,  as  expected. 

4.3  Ecosystem  simulations 

The  ROMS-EcoSim  model  with  its  baseline  initial  conditions  and  its  analytic  light  (AL)  model 
was  run  in  its  standard  mode  for  a  10-year  simulation  to  generate  the  baseline  ecosystem  prediction 
for  comparison  with  various  runs  using  the  EcoLight  (EL)  irradiance  calculations.  In  this  AL  run, 
the  irradiances  were  computed  at  every  grid  point  (denoted  in  the  figure  as  EGP),  every  ROMS  time 
step  (ETS),  every  wavelength  (5  nm),  and  down  to  the  Fa  =  0.001  (0.1%)  light  level.  Figure  3  shows 
the  resulting  surface  total  chlorophyll  values  over  the  10  years  (the  red  curve).  This  run  required 
5.28  hours  on  an  Apple  iMac  with  a  2.16  GHz  processor  and  1  Gbyte  of  RAM.,  or  about  32  minutes 
per  simulation  year. 

The  analytic  light  model  in  ROMS-EcoSim  was  then  replace  by  calls  to  EcoLight,  which  was 
run  in  various  modes  for  comparison  with  the  AL  baseline  run.  When  EcoLight  is  called  at  every 
grid  point,  every  time  step,  5  nm,  and  the  RTE  is  solved  down  to  the  F0  =  0.001  light  level,  the  run 
times  are  prohibitively  long  (estimated  to  be  400  hours  for  a  full  10-year  simulation).  However,  the 
EcoLight  chlorophyll  predictions  are  essentially  the  same  (to  within  1%  at  all  times  over  simulations 
of  a  year)  when  EcoLight  is  called  at  only  1  grid  point  (1GP).  The  irradiances  at  the  other  grid  points 
are  then  obtained  by  scaling  the  E0(z,e)  values  at  the  computed  grid  point  by  the  incident  sky 
irradiance  at  the  other  grid  points.  This  reduces  the  run  time  for  a  10-year  simulation  to  85.7  hours. 
This  EL,  1GP,  ETS,  5  nm,  Fa  =  0.001  simulation  is  taken  to  be  the  baseline  EcoLight  run  for 
comparison  with  the  10-year  AL  run  and  with  other  EL  runs.  The  resulting  surface  chlorophyll 
values  are  shown  in  green  in  Fig.  3. 

After  an  initial  period  of  numerical  and  biological  adjustment  to  the  initial  conditions,  both  light 
models  show  a  long-term  trend  of  slowly  decreasing  chlorophyll  values  from  one  year  to  the  next, 
with  annual-average  chlorophyll  values  (dotted  lines  in  Fig.  3)  decreasing  by  about  6%  a  year.  This 
is  likely  due  to  inadequate  vertical  mixing  to  replenish  the  nutrients  from  deep  water  below  the 
maximum  depth  of  the  simulation,  which  was  202  m.  It  should  be  remembered  that  these 
simulations  are  for  a  small  spatial  grid  that  cannot  reproduce  lateral  advection  as  occurs  in  a  true  3D 
ecosystem,  and  that  we  are  using  a  specific  parameter  setup  for  the  phytoplankton  populations. 
While  undesirable,  this  trend  does  not  affect  our  comparison  of  the  differences  in  ecosystem 
behavior  due  to  different  light  calculations,  with  all  else  being  the  same. 

Even  though  the  differences  are  not  great,  we  believe  that  this  baseline  EL  simulation  better 
represents  the  evolution  of  the  modeled  ecosystem  because  the  EcoLight  irradiances  are  accurately 
computed  throughout  the  water  column  at  all  wavelengths,  rather  than  being  obtained  from  the 
overly  simple  approximate  analytical  model  of  the  AL  simulation.  In  particular,  it  must  be 
remembered  that  the  phytoplankton  functional  groups  in  EcoSim  each  have  their  own  photosynthetic 
action  spectra,  so  that  they  utilize  different  wavelengths  in  different  ways,  and  thus  respond 
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differently  to  the  spectral  composition  of  the  scalar  irradiance  E0( z,e).  EcoLight  does  a  better  job 
of  modeling  E0( z,e)  near  400  nm,  beyond  600  nm,  and  at  depth  (Fig.  7,  discussed  below).  These 
differences  may  have  important  ecological  impacts  for  some  functional  groups. 

The  baseline  EL  simulation  run  time  is  still  far  too  long  for  routine  ecosystem  simulations. 
However,  it  is  not  necessary  to  recompute  the  irradiances  at  every  ROMS  time  step.  In  the  next  run, 
EL  was  called  at  1  grid  point  and  the  in-water  irradiance  was  recomputed  only  every  4  hours  (4HR), 
rather  than  at  every  ROMS  time  step.  The  most  recently  computed  irradiance  was  then  re-scaled  by 
the  input  sky  irradiance  at  the  current  time,  and  the  resulting  irradiance  profile  was  used.  The 
irradiances  were  computed  at  25  nm  resolution  and  interpolated  to  the  5  nm  resolution  required  by 
EcoSim.  Finally,  the  RTE  being  solved  down  to  only  the  15%  (F0  =  0.15)  light  level  and 
extrapolated  to  the  0.1%  level  needed  by  EcoSim  using  Eq.  (3).  With  these  optimizations,  the  total 
EL  run  time  was  reduced  to  6.88  hours.  The  resulting  ecosystem  predictions  remain  the  same  as  for 
the  baseline  EL  simulation  to  within  5%  at  all  times,  and  are  usually  closer.  This  is  an  increase  of 
only  30%  in  total  run  time  over  the  AL  simulation,  but  preserves  the  ecosystem  behavior  of  the  high- 
resolution  EL  simulation.  This  is  a  small  increase  in  total  computation  time  compared  to  the 
advantages  of  using  an  accurate  light  model  that  gives  a  wide  range  of  output  needed  for 
ecosystem  validation  or  for  initialization  and  data  assimilation. 

Figure  4  shows  year  ten  of  the  simulation.  The  AL  and  EL  baseline  runs  are  the  same  as  seen 
in  year  10  of  Fig.  3.  The  blue  curve  shows  the  optimized  EL  simulation  just  described:  1GP,  4HR, 
25  nm,  Fa  -  0.15. 

Another  run  was  made  with  EcoLight  further  optimized  to  1 GP,  6HR,  50  nm,  and  Fo  =  0.5 .  The 
resulting  run  time  decreased  to  6.73  hours.  However,  this  optimization  proved  to  be  too  coarse  of 
a  computational  scheme  in  that  the  resulting  ecosystem  behavior  differed  significantly  from  the  EL 
base.  The  surface  chlorophyll  values  differed  by  as  much  as  32%  (with  greater  differences  at  depth), 
and  the  time  of  the  spring  bloom  was  approximately  a  month  late.  These  results  are  thus  deemed  to 
be  unacceptable.  This  curve  is  shown  in  violet  in  Fig.  4. 

Figure  4  reveals  surprisingly  small  differences  between  the  long-term  AL  and  EL  annual  cycles 
(omitting  the  over-optimized  violet  EL  curve  from  further  consideration).  The  maximum  daily  near¬ 
surface  chlorophyll  difference  is  about  23%,  with  EL  being  greater  than  AL.  The  annual  averages 
are  2.65  mg  m  1  for  EL  and  2.42  mg  m  1  for  AL,  or  about  a  10%  difference.  The  times  of  the  spring 
bloom  are  about  the  same  in  both  models.  We  view  this  good  agreement  as  a  welcome  test  of  the 
accuracy  of  the  EcoSim  analytic  light  model  for  Case  1  water  and  as  a  model-model  validation  of 
the  correct  implementation  of  both  codes  for  computing  in- water  irradiances. 

Similar  results  are  found  at  greater  depths,  although  the  differences  can  be  somewhat  greater. 
Figure  5  shows  the  year  10  total  chlorophyll  values  at  depth  15.7  m.  At  this  depth,  the  maximum 
difference  between  AL  and  EL  is  36%,  although  the  annual  averages  are  within  4%  of  the  same. 
Figure  6  shows  the  depth  profiles  of  total  chlorophyll  at  day  I8O08  of  year  10,  when  the  water 
column  was  well  stratified  in  mid-summer. 
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The  run  parameters  and  the  run  times  for  these  simulations  are  summarized  in  Table  1.  Table 
2  shows  the  percent  of  total  run  time  used  for  the  physical  (ROMS)  vs.  bio-optical  (EcoSim  with  or 
without  EcoLight)  components  of  the  ecosystem  numerical  model. 


model 

time  resolution 

grid 

resolution 

F0 

wavelength 

resolution 

run  time 
in  hours 

Analytic 

9  min  (every 

ROMS  time  step) 

every  grid 
point 

0.001 

5  nm 

5.28 

EcoLight 
(run  1  year) 

9  min 

every  grid 
point 

0.001 

5  nm 

est.  400 

EcoLight 

9  min 

1  grid  point 

0.001 

5  nm 

87.5 

EcoLight 

every  4  hours 

1  grid  point 

0.15 

25  nm 

6.88 

EcoLight 

every  6  hours 

1  grid  point 

0.5 

50  nm 

6.73 

Table  1.  Options  used  in  comparison  runs.  The  mn  times  are  on  a  2.1  GHz  Mac  for  10-year 
simulations.  The  last  EcoLight  option  gave  unsatisfactory  results. 


simulation 

physical  module 

bio-optical  module 

AL  base  run 

60.2 

39.8 

EL  base  run 

3.7 

96.3 

EL  optimized  run  A 

46.2 

53.8 

Table  2.  Percent  of  total  run  time  for  10  year  simulations  used  by  the  physical  and  bio-optical 
modules  of  the  coupled  ecosystem  model. 


4.4  Irradiance  comparisons 

We  next  consider  how  the  analytic  and  numerical  light  models  compare  in  their  computed 
irradiance  values.  The  answer  can  be  seen  by  comparison  of  the  spectral  plane  irradiances  £d(z,e) 
at  time  0  of  the  simulation,  when  the  water  column  inherent  optical  properties  are  identical  in  all 
models.  Figure  7  shows  the  initial  Ed(z,e)  computed  by  the  AL  and  the  two  EL  models  used  in  Figs. 
4-6.  The  irradiances  for  these  three  models  are  close  throughout  the  euphotic  zone,  in  the  green  part 
of  the  spectrum.  However,  the  analytic -model  irradiances  are  less  the  EcoLight’ s  irradiances  near 
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400  nm  and  at  red  wavelengths,  with  the  differences  becoming  orders  of  magnitude  at  depth.  The 
two  EcoLight  versions  are  fairly  close  to  each  other  at  all  depths  and  wavelengths,  although  some 
differences  due  to  wavelength  interpolation  and  depth  extrapolation  are  apparent,  as  expected. 

The  reason  for  the  differences  in  the  analytic  and  EcoLight  spectral  scalar  irradiances  can  be 
traced  to  the  differences  in  the  mean  cosine  i  d  and  diffuse  attenuation  Kd  for  the  three  models.  A 
comparison  at  only  one  wavelength,  442  nm,  will  suffice.  Figure  8  shows  the  i  d(442)  depth  profiles 
at  time  0.  The  analytical  model  for  i  d(442)  is  much  different  than  the  mean  cosines  computed  by 
EcoLight.  In  particular,  the  analytical  i  d(442)  approaches  0.5  at  large  depths,  implying  an  isotropic 
downwelling  plane  irradiance,  which  is  not  correct  for  any  ocean  waters.  Even  at  asymptotic  depths 
in  homogenous  water,  the  radiance  distribution  is  far  from  isotropic  and  typical  i  d(442)  values  are 
around  0.7  to  0.8  at  great  depth,  as  seen  in  the  figure  for  EL  runA.  Note  that  EL  runB  solved  the 
RTE  only  to  the  15%  light  level,  which  was  about  10  m  at  this  time  and  wavelength.  Forcing  i  d  = 
0.5  at  great  depth  is  an  easily  correctable  weakness  of  the  default  analytic  mean  cosine  model  used 
in  EcoSim.  Note  that  in  the  full  EcoLight  calculation,  the  mean  cosine  continues  to  increase  slowly 
with  increasing  depth. 

In  the  analytic  irradiance  model,  the  mean  cosine  is  used  to  determine  the  Kd  depth  profile  as  seen 
in  Eq.  (2).  In  the  present  simulation,  i  d(442)  is  too  small  and  Kd(AA2)  is  consequently  too  large,  as 
seen  in  Fig.  9.  The  too-large  Kd  then  makes  Ed,  and  consequently  E0,  too  small,  as  seen  in  Fig.  7. 
These  analytic  vs.  EcoLight  differences  in  E0(z,e )  at  time  0  then  lead  to  slightly  different  biological 
concentrations  at  the  second  time  step.  These  differences  then  lead  greater  differences  in  ecosystem 
chlorophyll  values  as  time  goes  on. 


5.  Conclusions 


Figures  3-9  and  Tables  1  and  2  are  sufficient  to  establish  the  primary  results  of  this  study.  In 
spite  of  this  relatively  good  agreement  for  this  particular  idealized  simulation  of  deep  Case  1  water, 
we  argue  for  the  future  use  of  EcoLight  for  the  following  reasons: 

•  EcoLight  does  a  much  better  job  of  computing  irradiances  near  400  and  beyond  600  nm, 
and  at  great  depths.  These  differences  can  have  significant  effects  on  how  particular 
phytoplankton  functional  groups  evolve  in  time  and  at  depth. 

•  EcoLight  can  model  any  water  body,  including  Case  2  water  and  optically  shallow  waters 

for  which  bottom  reflectance  can  substantially  increase  the  irradiance  available  for 
photosynthesis  and  water  heating. 

•  EcoLight  provides  the  outputs  necessary  for  ecosystem  validation  using,  for  example, 
satellite-derived  water- leaving  radiances  or  remote-sensing  reflectances,  without  the  necessity 
of  using  a  chlorophyll  algorithm  to  convert  satellite  data  into  an  ecosystem  variable. 

•  EcoLight  provides  a  means  for  improved  ecosystem  model  initialization  and  data 
assimilation  directly  in  terms  of  the  available  data  (IOPs,  irradiances,  remote  sensing 
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reflectances),  again  without  an  intermediate  step  involving  a  data-to-chlorophyll  conversion 

algorithm. 

•  EcoLight  costs  at  most  30%  extra  in  the  total  simulation  run  time.  This  is  a  small 

computational  price  to  pay  for  its  advantages. 

Although  simple  analytic  irradiance  models  do  run  extremely  fast,  there  is  little  justification  for 
their  use  in  light  of  their  potential  inaccuracies  and  inconsistencies,  even  for  Case  1  waters. 
Moreover,  some  analytic  models  (although  not  EcoSim)  parameterize  the  water  inherent  optical 
properties  in  terms  of  a  single  quantity — the  total  chlorophyll  concentration — which  greatly 
oversimplifies  the  complex  relation  between  light  propagation  and  the  scattering  and  absorption 
properties  of  various  ocean  constituents.  Sophisticated  ecosystem  models  such  as  EcoSim  therefore 
track  several  functional  groups  of  particles,  each  of  which  has  its  own  absorption  and  scattering 
properties.  EcoSim  determines  its  total  inherent  optical  properties  as  sums  of  the  contributions  by 
any  number  and  type  of  components  and  is  thus  able  to  model  changes  in  the  light  field  induced  by 
changes  in  the  concentration  or  optical  properties  of  whatever  functional  groups  are  included  in  the 
biological  model.  Such  connections  between  ecosystem  components  and  the  light  field  cannot  be 
simulated  by  simple  analytic  models. 

Computational  stability  criteria  for  ecosystem  physical  models  may  require  time  steps  as  small 
as  one  minute.  However,  as  shown  here,  it  is  not  necessary  to  update  the  light  field  at  each  time  step. 
It  is  adequate  to  run  Ecolight  at  longer  intervals  and  then  interpolate  to  the  desired  time.  The  same 
idea  applies  to  running  Ecolight  on  a  coarser  spatial  grid  than  the  ecosystem  model  uses  and  then 
interpolate  to  obtain  the  required  spatial  resolution.  In  practice,  in  fully  3D  simulations,  it  would 
probably  be  possible  to  dynamically  determine  when  to  call  EcoLight.  For  example,  as  the 
calculations  proceed  to  a  new  grid  point,  the  IOPs  at  that  grid  point  could  be  compared  with  those 
at  the  last  grid  point  where  EcoLight  was  called  (the  reference  point).  If  the  IOPs  at  the  new  grid 
point  differ  by  less  than  some  amount  from  those  at  the  reference  point,  then  the  previously 
computed  irradiances  would  be  rescaled,  as  was  done  with  the  1GP  calculations  studied  here,  and 
applied  at  the  new  grid  point.  If  the  IOPs  at  the  new  grid  point  differ  by  more  than  some  amount 
from  those  at  the  reference  point,  then  EcoLight  would  be  called  for  a  cle  novo  calculation  of  the 
irradiance,  and  the  current  grid  point  would  become  the  new  reference  point.  A  dynamic 
determination  of  when  to  call  EcoLight  would  allow  for  frequent  calls  near  fronts  or  in  other 
situations  where  the  IOPs  or  external  properties  are  changing  rapidly,  and  for  few  calls  in  fairly 
homogenous  and  stable  ocean  regions.  When  employed  is  this  manner,  Ecolight  should  be 
applicable  to  three-dimensional  ecosystem  models  requiring  spectral  irradiances  at  many 
wavelengths  and  on  a  fine  spatial  and  temporal  grid,  but  at  little  more  computational  expense  than 
the  default  analytical  light  model  being  applied  at  every  time  step  and  grid  point,  as  is  done  now. 
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6.  Scientific  Impact 

Light  is  an  important  factor  in  determining  primary  production  and  mixed-layer  dynamics,  so 
EcoLight’s  accurately  computed  spectral  scalar  irradiances  will  significantly  improve  ecosystem 
model  predictions,  compared  to  the  use  of  approximate  analytical  models  for  spectral  irradiance  or 
PAR.  Moreover,  predictions  of  primary  production  in  shallow  waters  and  for  seagrass  beds,  for 
which  analytic  chlorophyll-based  models  can  be  off  by  a  factor  of  ten  or  more  near  highly  reflecting 
substrates,  would  be  well-served  by  the  use  of  Ecolight. 

In  addition,  EcoLight  computes  the  water-leaving  radiance  and  nadir-viewing  remote  sensing 
reflectance  at  no  extra  cost.  This  allows  for  direct  ecosystem  model  validation  using  satellite-derived 
remote- sensing  reflectances,  without  the  use  of  a  chlorophyll  algorithm  to  convert  the  satellite 
remote-sensing  reflectance  to  a  chlorophyll  value  for  comparison  withe  the  predicted  chlorophyll 
value.  Such  chlorophyll  algorithms  are  inaccurate  even  in  open-ocean,  Case  1  waters,  and  they  can 
be  in  error  by  over  an  order  of  magnitude  in  coastal,  Case  2  water.  They  are  not  applicable  at  all  in 
shallow  water  for  which  bottom  reflectance  influences  the  water- leaving  radiance. 

Finally,  EcoLight  computes  other  standard  optical  quantities  of  interest  such  as  various  in-water 
irradiances,  reflectances,  and  diffuse  attenuation  functions.  These  quantities  can  also  be  of  use  for 
ecosystem  validation  and,  importantly,  for  model  initialization  and  for  data  assimilation  when  data 
on  these  variables  are  available  from  observation  or  ancillary  models. 


7.  Future  Work 

We  have  properly  incorporated  EcoLight  as  a  subroutine  into  EcoSim  in  the  sense  that  EcoLight 
is  getting  the  correct  IOPs  and  other  information  from  ROMS  and  EcoSim,  solving  the  RTE 
correctly,  and  returning  accurate  irradiances  to  EcoSim  for  use  there.  However,  we  have  made  no 
attempt  to  rewrite  the  EcoLight  code  to  bring  it  up  to  ROMS-EcoSim  standards,  which  would  be 
necessary  for  efficient  simulations  on  large  spatial  scales  over  years  of  time.  The  EcoLight  code, 
which  is  called  as  a  subroutine  from  EcoSim,  remains  mostly  Fortran  77  and  contains  many  features 
that  are  not  really  needed  for  ecosystem  simulations.  We  speculate  that  a  thorough  rewriting  of 
EcoLight  in  Fortran  95  and  removal  of  unneeded  features  might  gain  another  factor  of  two  in  run 
times.  Use  of  a  faster  ODE  solver  for  depth  integration  of  the  RTE,  which  is  where  most  of 
EcoLight’s  run  time  is  used,  should  also  be  investigated. 
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simulation  year  1;  day  0.0 


Fig.  1.  Initial  depth  profiles  of  the  chlorophyll  concentrations  of  the  four  phytoplankton  functional 
groups,  and  of  the  total  chlorophyll. 


simulation  year  1 


Jan  Feb  Mar  Apr  May  Jun  Jul  Aug  Sep  Oct  Nov  Dec 
time  [calendar  month] 


Fig.  2.  Time  series  during  year  1  of  total  chlorophyll  at  1  m  depth  for  the  analytic  light  model  and 
3  sets  of  initial  conditions.  The  base  conditions  (red  curve)  at  time  0  are  those  seen  in  Fig.  1. 
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Fig.  3.  Ten-year  simulations  for  the  baseline  analytic  (AL,  red)  and  EcoLight  (EL,  green) 
simulations.  The  predicted  total  chlorophyll  values  are  plotted  at  local  noon  of  each  day.  The  annual 
average  Chi  values  are  shown  by  the  horizontal  dotted  lines  within  each  year. 
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Fig.  4.  Total  chlorophyll  values  at  1  m  depth  during  year  10  of  the  simulation.  The  irradiance 
calculation  models  are  identified  on  the  inset  labels. 
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simulation  year  10 
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Fig.  5.  Chlorophyll  concentrations  during  year  10  at  depth  of  15.7  m  for  various  simulations. 


simulation  year  10;  day180.0 


Fig.  6.  Chlorophyll  concentrations  at  as  a  function  of  depth  at  day  180  of  year  10. 
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wavelength  [nm] 


Fig.  7.  The  spectral  downwelling  plane  irradiances  Ed  at  time  0  at  selected  depths. 


Fig.  8.  The  mean  cosine  for  downwelling  radiance  at  442  nm,  i  d(442),  at  time  0.  Run  A  solved  the 
RTE  down  to  only  10  m  at  this  wavelength. 
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Kd(442)  [1/m] 


Fig.  9.  The  diffuse  attenuation  coefficient  for  downwelling  plane  irradiance  at  442  nm,  Kd{ 442),  at 
time  0.  The  Run  A  values  are  from  the  extrapolation  to  depth  by  Eq.  (3),  rather  than  from  solving 


the  RTE. 
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